#include <iterator>   // for stream_iterators, inserter
#include <iostream>
#include <fstream>
#include <string>
#include <boost/tr1/unordered_map.hpp>
#include <utility>
#include <algorithm>
#include <functional>
#include <limits>
#include <progress_indicator.h>
#include <boost/shared_ptr.hpp>
#include <boost/tr1/tuple.hpp>
#include <boost/range/algorithm.hpp>
#include <boost/format.hpp>
#include <boost/range/adaptor/map.hpp>
#include <boost/range/algorithm.hpp>
#include <boost/format.hpp>

#include <boost/filesystem.hpp>
#include <boost/filesystem/fstream.hpp>

using std::string;
using std::vector;
using std::ofstream;
using std::tr1::get;
using std::tr1::unordered_map;


#include "test_intervals_per_contig.h"
#include "parse_genomic_intervals.h"


// N.B. uses open-close coordinates

namespace test_intervals_per_contig
{
	namespace fs = boost::filesystem;

//________________________________________________________________________________________
//
//  get_next_interval
// 
//________________________________________________________________________________________
void t_test_intervals::get_next_interval(unsigned test_pos)
{
	//
	//	check not past the end 
	// 
	typedef vector<t_interval>::iterator i_intervals;
	using boost::range::lower_bound;
	i_intervals i_curr_interval = lower_bound(curr_contig_intervals->second, t_interval(test_pos, test_pos));
	if (i_curr_interval == curr_contig_intervals->second.end())
		curr_interval.beg = std::numeric_limits<unsigned>::max();
	else
		curr_interval = *i_curr_interval;
}

//________________________________________________________________________________________
//
//  get_flanking_intervals
// 
//________________________________________________________________________________________
void t_test_intervals::get_flanking_intervals(unsigned test_pos)
{
	//
	//	check not past the end 
	// 
	typedef vector<t_interval>::iterator i_intervals;
	using boost::range::lower_bound;
	i_next_interval = lower_bound(curr_contig_intervals->second, t_interval(test_pos, test_pos));

	if (i_next_interval ==)
	{
	}
	// 
	// entirely contained:
	// both flanking intervals point to the same thing
	// 
	if (test_pos >= curr_interval.beg && test_pos < curr_interval.end)
	{
		i_prev_interval = i_next_interval;
		return ;
	}

	// 
	// before next
	// 
	if (test_pos < curr_interval.beg)
		return false;


	if (test_pos < curr_interval.end)
		return true;

	if (i_next_interval == curr_contig_intervals->second.end())
	{
		curr_contig_intervals->second.end()
	}
	i_prev_interval = 
	if (i_curr_interval == curr_contig_intervals->second.end())
		curr_interval.beg = std::numeric_limits<unsigned>::max();
	else
		curr_interval = *i_curr_interval;
}


//________________________________________________________________________________________
//
//  c'tor
// 
//________________________________________________________________________________________
t_test_intervals::t_test_intervals(t_intervals_by_contig& intervals_by_contig_)
		:intervals_by_contig(intervals_by_contig_), curr_interval(0,0)
{
	//	
	//  Sort intervals
	// 
	using boost::range::sort;
	for (i_intervals_by_contig ii = intervals_by_contig.begin(), end = intervals_by_contig.end(); ii != end; ++ii)
		sort((*ii).second);
}

//________________________________________________________________________________________
//
//  contains
// 
//________________________________________________________________________________________
bool t_test_intervals::contains(const string& test_contig, unsigned test_pos)
{
	//
	//	new contig 
	// 
	if (curr_contig != test_contig)
	{	
		curr_contig = test_contig;
		curr_contig_intervals = intervals_by_contig.find(curr_contig);

		// 
		//  Missing contig
		//  	Set interval beginning to some really large value
		// 
		if (curr_contig_intervals == intervals_by_contig.end())
		{
			//cout << curr_contig << " not found\n";
			curr_interval.beg = std::numeric_limits<unsigned>::max();
			return false;
		}
		
		//
		//	Get next interval 
		// 
		else
			get_next_interval(test_pos);
	}
	//
	//  Is it within currect interval?
	// 
	if (test_pos < curr_interval.beg)
		return false;
	if (test_pos < curr_interval.end)
		return true;


	// 
	// 	Beyond this interval: get next
	// 
	get_next_interval(test_pos);

	//
	//  Is it before or within currect interval?
	//  	Cannot be after: otherwise last interval, and curr_interval.beg == MAX_UINT
	// 
	return test_pos >= curr_interval.beg;
}

//________________________________________________________________________________________
//
//  contains
// 
//________________________________________________________________________________________
// returns distance and interval_id
std::pair<unsigned, unsigned> t_test_intervals::closest(const string& test_contig, unsigned test_pos)
{
	//
	//	new contig 
	// 
	if (curr_contig != test_contig)
	{	
		curr_contig = test_contig;
		curr_contig_intervals = intervals_by_contig.find(curr_contig);

		// 
		//  Missing contig
		//  	Set interval beginning to some really large value
		// 
		if (curr_contig_intervals == intervals_by_contig.end())
		{
			//cout << curr_contig << " not found\n";
			unsigned max_uint = std::numeric_limits<unsigned>::max();
			curr_interval.beg = max_uint;
			return std::make_pair<unsigned, unsigned(max_uint, max_uint);
		}
		
		//
		//	Get next interval 
		// 
		else
			get_next_interval(test_pos);
	}
	//
	//  Is it within currect interval?
	// 
	if (test_pos < curr_interval.beg)
		return false;
	if (test_pos < curr_interval.end)
		return true;


	// 
	// 	Beyond this interval: get next
	// 
	get_next_interval(test_pos);

	//
	//  Is it before or within currect interval?
	//  	Cannot be after: otherwise last interval, and curr_interval.beg == MAX_UINT
	// 
	return test_pos >= curr_interval.beg;
}


void read_intervals(t_intervals_by_contig& intervals_by_contig, fs::ifstream& f_regions_to_retain)
{
	unsigned cnt_line = 0;

	sp_genomic_interval_grammar genomic_interval_grammar = create_genomic_interval_grammar(cnt_line);

	string str;
	// ignore first line
	while (std::getline(f_regions_to_retain, str))
	{
		cnt_line ++;
		if (!str.length() || str[0] == '#')
			continue;

		t_genomic_interval interval;
		if (!parse_genomic_interval(genomic_interval_grammar, str, interval))
			break; 

		// convert to open / closed by adding 1
		intervals_by_contig[interval.contig].push_back(t_interval(interval.beg, interval.end + 1));
	}
}



}
